home *** CD-ROM | disk | FTP | other *** search
/ Chip 2007 January, February, March & April / Chip-Cover-CD-2007-02.iso / Pakiet bezpieczenstwa / mini Pentoo LiveCD 2006.1 / mpentoo-2006.1.iso / livecd.squashfs / usr / lib / python2.4 / site-packages / Numeric / Matrix.py < prev    next >
Encoding:
Python Source  |  2006-01-20  |  6.0 KB  |  213 lines
  1. import string
  2. import types
  3.  
  4. import Numeric
  5. import LinearAlgebra
  6. from UserArray import UserArray, asarray
  7. from Numeric import dot, identity, multiply
  8. from MLab import squeeze
  9.  
  10. # make translation table
  11. _table = [None]*256
  12. for k in range(256):
  13.     _table[k] = chr(k)
  14. _table = ''.join(_table)
  15.  
  16. _numchars = string.digits + ".-+jeEL"
  17. _todelete = []
  18. for k in _table:
  19.     if k not in _numchars:
  20.         _todelete.append(k)
  21. _todelete = ''.join(_todelete)
  22.  
  23. def _eval(astr):
  24.     return eval(astr.translate(_table,_todelete))
  25.  
  26. def _convert_from_string(data):
  27.     data.find
  28.     rows = data.split(';')
  29.     newdata = []
  30.     count = 0
  31.     for row in rows:
  32.         trow = row.split(',')
  33.         newrow = []
  34.         for col in trow:
  35.             temp = col.split()
  36.             newrow.extend(map(_eval,temp))
  37.         if count == 0:
  38.             Ncols = len(newrow)
  39.         elif len(newrow) != Ncols:
  40.             raise ValueError, "Rows not the same size."
  41.         count += 1
  42.         newdata.append(newrow)
  43.     return newdata
  44.  
  45.  
  46. _lkup = {'0':'000',
  47.          '1':'001',
  48.          '2':'010',
  49.          '3':'011',
  50.          '4':'100',
  51.          '5':'101',
  52.          '6':'110',
  53.          '7':'111'}
  54.  
  55. def _binary(num):
  56.     ostr = oct(num)
  57.     bin = ''
  58.     for ch in ostr[1:]:
  59.         bin += _lkup[ch]
  60.     ind = 0
  61.     while bin[ind] == '0':
  62.         ind += 1
  63.     return bin[ind:]
  64.         
  65.  
  66. class Matrix(UserArray):
  67.     def __init__(self, data, typecode=None, copy=1, savespace=0):
  68.         if type(data) is types.StringType:
  69.             data = _convert_from_string(data)
  70.         UserArray.__init__(self,data,typecode, copy, savespace)
  71.         if len(self.array.shape) == 1:
  72.             self.array.shape = (1,self.array.shape[0])
  73.         self.shape = self.array.shape
  74.  
  75.     def __getitem__(self, index):        
  76.         out = self._rc(self.array[index])
  77.         try:
  78.             n = len(index)
  79.             if n > 1 and isinstance(index[0], types.SliceType) \
  80.                and isinstance(index[1], types.IntType):
  81.                 sh = out.array.shape
  82.                 out.array.shape = (sh[1],sh[0])
  83.                 out.shape = out.array.shape
  84.         except TypeError:  # Index not a sequence
  85.             pass
  86.         return out
  87.  
  88.     def __mul__(self, other):
  89.         aother = asarray(other)
  90.         if len(aother.shape) == 0:
  91.             return self._rc(self.array*aother)
  92.         else:
  93.             return self._rc(dot(self.array, aother))
  94.  
  95.     def __rmul__(self, other):
  96.         aother = asarray(other)
  97.         if len(aother.shape) == 0:
  98.             return self._rc(aother*self.array)
  99.         else:
  100.             return self._rc(dot(aother, self.array))
  101.  
  102.     def __imul__(self,other):
  103.         aother = asarray(other)
  104.         self.array = dot(self.array, aother)
  105.         return self
  106.  
  107.     def __pow__(self, other):
  108.         shape = self.array.shape
  109.         if len(shape)!=2 or shape[0]!=shape[1]:
  110.             raise TypeError, "matrix is not square"
  111.         if type(other) in (type(1), type(1L)):
  112.             if other==0:
  113.                 return Matrix(identity(shape[0]))
  114.             if other<0:
  115.                 result=Matrix(LinearAlgebra.inverse(self.array))
  116.                 x=Matrix(result)
  117.                 other=-other
  118.             else:
  119.                 result=self
  120.                 x=result
  121.             if other <= 3:
  122.                 while(other>1):
  123.                     result=result*x
  124.                     other=other-1
  125.                 return result
  126.             # binary decomposition to reduce the number of Matrix
  127.             #  Multiplies for other > 3.
  128.             beta = _binary(other)            
  129.             t = len(beta)
  130.             Z,q = x.copy(),0
  131.             while beta[t-q-1] == '0':
  132.                 Z *= Z
  133.                 q += 1
  134.             result = Z.copy()
  135.             for k in range(q+1,t):
  136.                 Z *= Z
  137.                 if beta[t-k-1] == '1':
  138.                     result *= Z
  139.             return result
  140.         else:
  141.             raise TypeError, "exponent must be an integer"        
  142.  
  143.     def __rpow__(self, other):
  144.         raise TypeError, "x**y not implemented for matrices y"
  145.  
  146.     def __invert__(self):
  147.         return Matrix(Numeric.conjugate(self.array))
  148.     
  149.     def __setattr__(self,attr,value):
  150.         if attr=='shape':
  151.             if len(value) == 0:
  152.                 value = (1,1)
  153.             if len(value) == 1:
  154.                 value = (1,) + value 
  155.             if len(value) != 2:
  156.                 raise ValueError, "Can only reshape a Matrix as a 2-d array."
  157.             self.array.shape=value
  158.         self.__dict__[attr]=value
  159.  
  160.     def __getattr__(self, attr):
  161.         if attr == 'A':
  162.             return squeeze(self.array)
  163.         elif attr == 'T':
  164.             return Matrix(Numeric.transpose(self.array))
  165.         elif attr == 'H':
  166.             if len(self.array.shape) == 1:
  167.                 self.array.shape = (1,self.array.shape[0])
  168.             return Matrix(Numeric.conjugate(Numeric.transpose(self.array)))
  169.         elif attr == 'I':
  170.             return Matrix(LinearAlgebra.inverse(self.array))
  171.         elif attr == 'real':
  172.             return Matrix(self.array.real)
  173.         elif attr == 'imag':
  174.             return Matrix(self.array.imag)
  175.         elif attr == 'flat':
  176.             return Matrix(self.array.flat)
  177.         else:
  178.             raise AttributeError, attr + " not found."
  179.  
  180.     def __setitem__(self, index, value):
  181.         value = asarray(value, self._typecode)
  182.         self.array[index] = squeeze(value)
  183.     def __setslice__(self, i, j, value): self.array[i:j] = asarray(squeeze(value),self._typecode)
  184.  
  185.     def __float__(self):
  186.         return float(squeeze(self.array))
  187.     
  188.     def m(self,b):
  189.         return Matrix(self.array * asarray(b))
  190.     
  191.     def asarray(self,t=None):
  192.         if t is None:
  193.             return self.array
  194.         else:
  195.             return asarray(self.array,t)
  196.  
  197. if __name__ == '__main__':
  198.         from Numeric import *
  199.         m = Matrix( [[1,2,3],[11,12,13],[21,22,23]])
  200.         print m*m
  201.         print m.array*m.array
  202.         print transpose(m)
  203.         print m**-1
  204.         m = Matrix([[1,1],[1,0]])
  205.         print "Fibonacci numbers:",
  206.         for i in range(10):
  207.             mm=m**i
  208.             print mm[0][0],
  209.         print
  210.         
  211.             
  212.             
  213.